#!/usr/bin/python
"""
Usage: python auto_anno.py annofile.txt seqfile.fasta
Author:			Jimmy Saw
Date modified:	12-20-2011

Description:    This program parses Illumina reads aligned against an assembled contig.

Usage:          python dissertation_IlluminaCoverage.py <coordinate file> <MAXCOV>
                MAXCOV = maximum number of reads to cover a given span.

Example:        python dissertation_IlluminaCoverage.py contig.b1.coords 100

Note:           Needs MUMMER output file generated like this:
nucmer contig.fasta B_trimmed.1.fasta -p contig.b1
show-coords -r -c -l -T contig.b1.delta > contig.b1.coords

"""

import sys
import re
from Bio import SeqIO

mummer_file = sys.argv[1]
mf = open(mummer_file, "rU")
mfl = mf.readlines()

MAXCOV = int(sys.argv[2]) #maximum number of Illumina reads to recruit

allreads = []

for line in mfl[4:]:
    strand = ""
    l = line.split('\t')
    ref_start = int(l[0])
    ref_stop = int(l[1])
    qry_start = int(l[2])
    qry_stop = int(l[3])
    if qry_stop > qry_start:
        strand = "+"
    else:
        strand = "-"
    per_id = float(l[6])
    qry_cov = float(l[10])
    qry_name = l[12].rstrip()
    allreads.append((ref_start, ref_stop, qry_start, qry_stop, strand, 
                   per_id, qry_cov, qry_name))

tmpstr = mfl[4].split('\t')
contig_size = int(tmpstr[7])

to_keep_list = []

i = 0
while i < contig_size:
    readlist = []
    for read in allreads:
        rstart = read[0]
        rstop = read[1]
        if 0+i < rstart < i+50:
            readlist.append(read)
    if len(readlist) > MAXCOV:
        for j in readlist[0:5]:
            to_keep_list.append(j)
    else:
        for j in readlist:
            to_keep_list.append(j)
    i += 80 #increment by 80 so that the next read overlaps with previous read by about 20bp

for i in to_keep_list:
#    print i[-1], i[0]
    print i[-1]

mf.close()

